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We apply the concept of multistep cascading to the problem of fourth-harmonic generation in 
a single quadratic crystal. We analyze a new model of parametric wave mixing and describe its 
stationary solutions for two- and three-color plane waves and spatial solitons. Some applications to 
the optical frequency division as well as the realization of the double-phase-matching processes in 
engineered QPM structures with phase reversal sequences are also discussed. 



Cascading effects in optical materials with quadratic 
(second-order or x^) nonlinear response provide an effi- 
cient way to lower the critical power of all-optical switch- 
ing devices [Q. The concept of multistep cascading ||] 
brings new ideas into this field, leading to the possibil- 
ity of an enhanced nonlinearity-induced phase shift and 
generation of multicolor parametric spatial solitons. In 
particular, multistep cascading can be achieved by two 
nearly phase matched second-order nonlinear processes, 
second-harmonic generation (SHG) and sum-frequency 
mixing (SFM), involving the third-harmonic wave 
In this Letter, we extend the concept of multistep cas- 
cading to nonlinear effects of the fourth order and the 
fourth-harmonic generation (FHG) in a single noncen- 
trosymmetric crystal. In particular, we analyze a new 
model of multistep cascading that involves the FHG pro- 
cess, and describe its stationary solutions for normal 
modes — plane waves and spatial solitons. Our study 
provides the first systematic analysis of the problem of 
FHG via a pure cascade process, observed experimen- 
tally more than 25 years ago || and later studied in a 
cascading limit only Bj. 

We consider the FHG via two second-order parametric 
processes: to + uj — 2lu and 2uj + 2uj — Alu, where to is the 
frequency of the fundamental wave. In the approxima- 
tion of slowly varying envelopes with the assumption of 
zero absorption of all interacting waves, we obtain 



+ ij 1 A*Se~ iAklZ , 

+ i lx A 2 e tAklZ + i l2 S*Te- tAk2Z , 

+ i l2 S 2 e lAk2Z , 
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where A, S, and T are the envelopes of the 
fundamental- frequency (FF), second- (SH) and fourth- 
harmonic (FH) waves respectively, 71,2 are propor- 
tional to the elements of the second-order suscep- 
tibility tensor, and Afei = 2ki — k 2 and Ak 2 — 
2k 2 — &4 are the corresponding wave-vector mis- 
match parameters. We introduce the normalized en- 
velopes (it, v, w) according to the following relations: 



A(x,z) = (16z d "/ 2 )~ 1 u(x/a 1 z/2z d )exp(-iAkiz/2), 
S(x,z) — (8z d 'y 2 )^ 1 v(x/a, z/2zd), and T(x.z) = 
(4zd'j 2 )~ 1 w(x/a, z/2zd) exp(iAk 2 z), where a is the char- 
acteristic beam width, and Zd — k\a 2 is the diffraction 
length of the FF component. In order to describe a fam- 
ily of nonlinear modes characterized by the 'propagation 
constant A, we look for solutions in the form u(x, z) — > 
XU(xy/\\\,z\X\)e lXz / 4 , v(x,z) -► XV{xy/\X\, z\\\)e iXz ^ 2 , 
and w(x,z) — > AW(x\/[A[, z|A|)e lAz , and obtain the nor- 
malized equations: 
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ai U + X U*V = 0, 
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\~'W+^U 2 = 0, 
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(1) 



aW + -V z = 0. 



Here s = sign(A) = ±1, \ = 71/(472) is a relative 
strength of two parametric processes, and the normal- 
ized mismatches are defined as a = 4 + [3/\ and a»i = 
1/4 + Pi/X, where [3 = 8Ak 2 z d and /3i = -Akiz d . 

First, we analyze the plane- wave solutions of Eq. (|l|) 
which do not depend on x. In this case, the total intensity 
I is conserved, and we present it in terms of the unsealed 
variables as I = I u +I v +I w , where I u — \u\ 2 /A, I v = \v\ 2 , 
and I w = 4|u;| 2 . Solutions {Uo, Vo, Wo}, which do not de- 
pend on z, are the so-called normal modes. The simplest 
one-component FH mode {0,0, Wo} exists at a — 0. It 
has a fixed phase velocity A = — /3/4 and an arbitrary 
amplitude, being unstable for I w > /? 2 /4 due to a para- 
metric decay instability. 

Two-mode solution {0, \f2~a, 1} describes a paramet- 
ric coupling between SH and FH waves, and it exists for 
a > 0, bifurcating at a — from the FH mode. Cou- 
pling of this two-mode plane wave to a FF wave can lead 

(cr) 



to its decay instability, provided |ai| < af 1 ' = xV2a 
To understand the physical meaning of this inequality, 
we note that the family of solutions characterized by 
the propagation constant A corresponds to a straight 
line in the (a, ai)-parameter space, see Figs. |l|(a,b). 
Moreover, all such lines include the point (4, 1/4) as the 
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asymptotic limit for |A| — > +00. This special point be- 
longs to the instability region if the relative strength of 
the FF-SH interaction exceeds a critical value, i.e. for 
x > X M = 1/(8^/2) ~ 0.088. However, for X < X (cr) this 
decay instability is suppressed for highly intense waves, 
owing to a strong coupling with the FH field. 
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FIG. 1. (a,b) Existence and stability of three-mode plane 
waves for (a) \ = x' cr) /4 and (b) x = 4x (cr) - Light shad- 
ing — stable, black — unstable, dark shading — oscillatory 
unstable, and blank — no solutions. Open circles mark ex- 
act phase- matching, dash-dotted lines correspond to the lower 
plots. (c,d) Intensity vs. A for (c) /3 = (3% = —0.5; (d) = — 1 
and /3i = 0.04. Thick dashed/solid vertical line — one-wave 
(FH) modes, thick curves — two-wave (SH + FH) modes, 
thin curves — three-wave modes. Solid — stable, dashed — 
unstable, and dotted — oscillatory unstable modes. Open cir- 
cles mark the bifurcation points. The legend beneath shows 
the stability of three-wave modes by using the same shadings 
as in the upper plots. 



Finally, a three-mode solution, 
V 2 /(2a), U = y / 2Va(T 



and < a.\ < 
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Vq = ai/x, W = 
exists for (i) a > 

and 
such 
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W Q )/x, 
(ii) a > and a\ < —ct( 
(hi) a < and a± > 0. In the limit |A| — * +00, 
three-wave modes are possible only for x > X 
the region (i), stability properties of the three- wave so- 
lutions are determined by a simple criterion: the modes 
are stable if dI/d\X\ > 0, and unstable, otherwise. For 
the parameter regions (ii) and (iii), oscillatory instabil- 
ities are possible as well. Existence and stability of all 
types of stationary plane- wave solutions of the model (|l|) 
are summarized in Figs. |l|(a-d). 

In general, the system (Q) is nonintegrable and its dy- 
namics are irregular. However, we find that in some cases 
a quasi-periodic energy exchange between the harmonics 
is possible. Figure ^(a) shows one such case, when the in- 
tensities of unstable two-wave and stable three- wave sta- 
tionary modes are close to each other, and an unstable 
two- wave mode periodically generates a FF component. 
Less regular dynamics are observed for other cases, such 
as for the generation of both SH and FH waves from an 



input FF wave [Fig. g(b)]. This example also illustrates 
the possibility of effective energy transfer to higher har- 
monics close to the double phase matching point. 





FIG. 2. Dynamics of plane waves: (a) Instability 
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a two- wave mode [corresponding to A = 2.1 in Fig. Mc)], 
(b) Generation of higher harmonics from a FF input [param- 
eters correspond to Fig. [j](d)]. FF, SH, and FH components 
are shown by dotted, solid, and dashed curves, respectively. 

Equations (0) may have a different physical meaning 
provided the normalized amplitude v stands for the mode 
of the fundamental frequency u>. Then, Eqs. (Q) de- 
scribe the optical frequency division by two (the field u) 
via parametric amplification and down-conversion (see, 
e -g-i @); provided both FF (v) and SH (w) fields are 
launched simultaneously at the input. Such a frequency 
division parametric process is in fact shown in Fig. ^(a), 
where this time the generated u wave corresponds to the 
frequency u>/2, and is shown as dotted. 

We now look for spatially localized solutions of 
Eqs. (Q), quadratic solitons. First, we note that two- wave 
solitons consistingof the SH and FH components can be 
approximated as : 

V (x) = V m sech p (x/p), Wo(x) - W m sech 2 (a;/p), (2) 
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where all parameters are functions of a only. Bright soli- 
tons either do not exist or are unstable being in reso- 
nance with linear waves [B] outside the parameter region 
a > and a± > 0, at s = +1. We find that in this re- 
gion, similar to plane-wave modes, three-wave solitons 
exist for < ct\ < a^, where the critical (cut-off) 



value ct^ corresponds to a bifurcation from the two- 
wave solution (^]). In order to find a^ CI \ we should solve 
the first equation of the system (Q) with the SH profile 
from Eq. (|^) (sec also ||). However, such a linear eigen- 
value problem has no exact analytical solution for arbi- 
trary p, and thus we introduce an approximation Vo(x) ~ 
Vq(x) = VmSech 2 (x / q) , requiring that the functions co- 
incide at the amplitude level V m /2, and then define the 
scaling parameter as q — p cosh~ 1 (2 1 / p )/cosh _1 (2 1/ ' 2 ). 
Such an approximation adequately describes an effective 
soliton waveguide, and thus should provide overall good 
accuracy (except for some limiting cases). After solving 
the eigenvalue problem with the potential xVq(x), we ob- 
tain an approximate expression for the bifurcation points: 
a[" } ~ + W mX q 2 - 1 - 2n) 2 /(4q 2 ), where n is the 
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order of the mode guided by the two-component paramet- 
ric soliton waveguide [0] . For a single- hump mode (n = 0) 
the behavior of this cut-off is very similar to that of the 
plane waves. Indeed, in the cascading limit (a 3> 1), we 
have V m ~ 2y / a, and ~ 2x\fa, which differs by y/2 
from the corresponding result for plane waves. The crit- 
ical value of x f° r one-hump solitons can also be found 
from the approximate solution, x ~ 0.132. We per- 
formed numerical simulations and found that the accu- 
racy of our approximation is of the order of (and usually 
better than) 1% in a wide range of parameters (\ > 1 0~ 2 
and a > ICT 2 ), see Fig. ||(a). 




FIG. 3. (a) Regions of existence and stability of 
three-mode parametric solitons [shading is the same as in 
Fig. |l|(a)]. Open diamonds — an analytical approximation, 
dark circle — exact phase-matching point. The dash-dotted 
line corresponds to the solutions at f3 = 2 and fii = —0.15, for 
which the power vs. A dependences are shown in (b): thick — 
two- wave (SH + FH), and thin — three- wave solitons; solid 
and dashed lines mark stable and unstable solutions, respec- 
tively. Open circle is the bifurcation point. (c,d) Development 
of a decay instability of a two-wave soliton corresponding to 
A = 1 in (b), and generation of a three-component soliton: 
(c) FF, SH, and FH peak intensities vs. distance shown by 
dotted, solid, and dashed curves, respectively; (d) evolution 
of the FF component. For all the plots \ — X /2. 

Quite remarkably, for both positive a and a\ the sta- 
bility properties of solitons [see Figs. ||(a,b)] and plane 
waves [see Figs. |f|(a,c)] look similar. Specifically, sta- 
bility of two- and three-component solitons is defined 
by the Vakhitov-Kolokolov criterion dP/dX > 0, where 
P = f_°° Idx is the soliton power, except for the region 

(cr) 

d\ < a\ where two-component solutions exhibit para- 
metric decay instability. An example of such an instabil- 
ity is presented in Figs. ||(c,d), where an unstable two- 
wave soliton generates a stable three-wave state. Such 
instability-induced dynamics are very different from that 
of plane waves where, instead, quasi-periodic energy ex- 
change is observed [see Fig. §(a)] . In the case of localized 
beams, diffraction leads to an effective power loss and 



convergence to a new (stable) state. 

Similar to other models of multistep cascading P,Q, 
Eqs. (|l|) possess various types of exact analytical solu- 
tions, which can be found at a = ot\ = 1 and x > 
1/V2; ax = a/4 (0 < a < 1) and X = 1/0/2) or 
X = [(3a)/(4 + 2a)] 1 / 2 . Details will be presented else- 
where. 

In order to observe experimentally the multistep 
cascading and multi-frequency parametric effects de- 
scribed above, we should satisfy the double-phase match- 
ing conditions. Using the conventional quasi-phase- 
matching (QPM) technique || for FHG via a pure cas- 
cade process in LiTa03, we find that there exists only 
one wavelength (2.45^tm), for which two parametric pro- 
cesses can be phase-matched simultaneously by the dif- 
ferent orders m of the QPM structure with the pe- 
riod Aq = 34/im. However, for the so-called phase- 
reversal QPM structures |l(| characterized by two pe- 
riods, the QPM period Aq and the modulation period 
A p h (A p h > Aq), double-phase matching is possible in 
a broad spectral range, provided the periods are se- 
lected to satisfy the conditions: Aq = 47r|Afci + Afc 2 | _1 , 
A p h = 47rm|Afci — A/c 2 p , where m is the grating order. 
Thus, the engineered QPM structures suggested in |l(]] 
are more efficient than, e.g., the Fibonacci superlattices, 
and they can be used to achieve double-phase match- 
ing and to support different types of multistep cascading 
processes. 

In conclusion, we have introduced a new model of 
multistep cascading that describes the fourth-harmonic 
generation via parametric wave mixing. We have ana- 
lyzed the existence and stability of the stationary solu- 
tions of this model for normal modes — plane waves and 
spatial solitons. We have also discussed the possibility 
of double-phase-matching in engineered QPM structures 
with phase-reversal sequences. 
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